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Abstract 



^ • In Part I of this paper we have presented a simple model capable of describing 

O . the localized failure of a massive structure. In this part, we discuss the identification 

of the model paiameters from two kinds of experiments: a uniaxial tensile test and a 
three-point bending test. The former is used only for illustration of material parame- 
ter response dependence, and we focus mostly upon the latter, discussing the inverse 
\^ ', optimization problem for which the specimen is subjected to a heterogeneous stress 

^ ■ field. 



(N 

^ 1 Introduction 

o 

Present work is the follow-up on Part I of this paper [[11, discussing a relatively simple 
I model capable of describing the behavior of a massive structure until the point of local- 

^ I ized failure. The model contains all the ingredients for taking into account both the diffuse 

damage mechanism, which leads to the appearance of microcracks, as well as the failure 
process characterized by the propagation of macrocracks. Perhaps the most important 
advantage of the proposed model is the fact that all its parameters have a clear physical 
interpretation and can be straightforwardly visualized in terms of the shape of a stress- 
strain diagram. In addition, influence of each parameter is dominant only for a specific, 
easily recognizable, stages of material behavior. This kind of a priori knowledge has a po- 
tential to greatly simplify the model calibration and will be systematically used throughout 
the paper. 

In Part II of this paper, we discuss the identification of the model parameters from 
experimental measurements made on a structural level. Generally speaking, the com- 
plexity of the identification procedure is determined by the choice of experimental setup. 
Solely from the identification point of view, the simplest experiment to execute is the 
uniaxial tensile test. In this case, the strain field stays mostly homogeneous during the 
whole procedure and the global response represented by the load-displacement diagram 
is very similar to the stress-strain curve for one material point; see Section [3] for more 
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details. The model parameters can be then directly determined from the shape of the 
load-displacement curve. Such a uniform loading is, however, very difficult if not impos- 
sible to impose in a laboratory test, especially for quasi-brittle materials. Therefore, other 
tests are often used in experimental practice. 

The three-point bending test, in particular, is considered to be much simpler to perform 
and its results are well-reproducible. Therefore, we focus on the identification procedure 
for the proposed model parameters directly from results of three-point bending test. Main 
difficulty is in this case imposed by heterogeneity of the stress and the strain fields, which 
is present since the very start of the experiment. The macro-scale measurements pro- 
vide the load-deflection curve that integrates data from different parts of the specimen 
experiencing different regimes of (in)elastic behavior. For that reason, the possibility of 
a simple determination of model parameters from load-deflection curve is lost and an 
advanced calibration procedure needs to be applied. 

To take advantage of the model specific structure, already mentioned above, the identi- 
fication procedure should be divided into three sequential stages discussed in detail in Sec- 
tionlH From the algorithmic point of view, the material calibration can then be understood 
as a sequential optimization problem. Such approach has two main advantages: first, solv- 
ing three simpler identification steps in a batch form is typically much more efficient then 
the full-scale problem; second, it allows to use only a subset of simulations for initial 
stages of identification process. 

A variety of techniques is available to identify material parameters via optimization 
methods, see e.g. [19. and reference therein]. The gradient-based methods are usually 
considered to be the most computationally efficient optimization algorithms available 
and as such have been successfully used in a variety of identification problems, e.g. 
[fT8l [TOl l20ll . For the current model, however, analytic determination of sensitivities is 
fairly difficult, mainly due the history dependency of the model as well as complex inter- 
action of individual parameters. The accuracy of numerical approximation to the 'exact' 
sensitivities, on the other hand, is driven by the choice of pseudo-time step used in numer- 
ical simulations. Clearly, to reduce the computational time, the pseudo-time step should 
be used as large as possible. Therefore, the response-based objective function will not be 
smooth and gradient-based methods are unlikely to be very successful. 

As an alternative, techniques of soft-computing can be employed for optimization 
of complex objective functions. For example, stochastic evolutionary algorithms have 
been successfully used for solution of identification problems on a level of material point 
[f5l l25ll or on a level of simple structures lITTl ITTll . For the current case, however, com- 
plexity of the optimization can be attributed rather to its non- smooth character than to the 
appearance of multiple optima; the family of problems where evolutionary algorithms are 
the most successful methods. This opens the way to more specialized tools, which deliver 
higher efficiency when compared to usually time-consuming evolutionary algorithms. 

The approach adopted in the present work is based on an adaptive smoothing of the ob- 
jective function by artificial neural networks (see, e.g., [|24ll26ll for alternative ANN-based 
solutions to identification problems). In particular, the approximated model is provided 
by the Radial Basis Function Network, described in Section 15.21 dynamically evolved 
by minima located by a real-encoded genetic algorithm, briefly reviewed in Section 15. 1[ 
The proposed sequential numerical strategy is systematically verified in Section [6] with at- 
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tention paid to a detailed assessment of the proposed stochastic algorithm reliability. Final 
remarks and conclusions can be found in Section |71 



2 A brief description of the identified model 

In the present section, we give a brief description of the model on which the identification 
procedure is based. For the readers interested in more details, the complete description 
of the model is given in Part I of this paper. 

As already mentioned, the proposed model is capable of taking into account two dif- 
ferent types of dissipation (e.g. see 1121 ): 

• a bulk dissipation induced by the appearance of uniformly distributed microcracks. 
This bulk dissipation is taken into account by the use of a classical continuum dam- 
age model; 

• a surface dissipation induced by the development of macrocracks responsible for the col- 
lapse of the structure. As presented in Part I of this paper, this phase is taken into 
account by the use of a strong discontinuity model. The surface dissipation is taken 
into account by the introduction of a traction/displacement jump relation. 

Therefore, two different models are involved in the constitutive description: the one 
associated with the bulk material and the one associated with the displacement disconti- 
nuity. Both are built on the same scheme considering the thermodynamics of continuous 
media and interfaces. 

The key points of the construction of each of the two models are summarized in Table 
[Hand Table II 



Helmholtz energy 




Damage function 


^{(T,q) = V(T:D':(T--j^{Gf-q) 


State equations 


(T = D^:e and q = —^Z{^) 


Evolution equations 




Dissipation 





Table 1 : Main ingredients of the continuum damage model 



For the discrete damage model, the isotropic softening law is chosen as: 

1 — exp 



(1) 



In Tables [U and |2l the variables y, y\ and 72 denote Lagrange multipliers induced 
by the use of the maximum dissipation principle. Moreover, u denotes the displacement 



3 



Helmholtz energy 


^(a,Q,l) = ia-0-^-a+g(|) 


Damage functions 


_ ^Ktr,,!) =tr,-n-(a-/-|) 
l2(tr„^) = |tr,-m|-(^.-f'?) 


State equations 


tr^=Q-i.a and ^ = -§ 


Evolution equations 


^-7itr;.„+r2|t,;.„,i ; 1-71 + J; 72 


Dissipation 


0^# = i|(a/-l|) 



Table 2: Main ingredients of the discrete damage model 



E 


e 


(25.0, 50.0) GPa 


V 


G 


(0.1,0.4) 




e 


(1.0,5.0) MPa 


K 


e 


(10.0, 10000.0) MPa 




e 


(O-/+0.1, 2o-/) MPa 




G 


(0.1 a/, lO.Oo-/) MPa/mm 



Table 3: Limits for the model parameters. 



jump on the surface of discontinuity. Finally, D and Q correspond to the damaged com- 
pliance of the continuum and discrete model, respectively. 

Note that in a three-point bending or a uniaxial tensile test, the simulated response 
is almost independent of the limit tangential traction Os- Therefore, its value was set 
to O.ldf. With such simplification, there are six independent material parameters to be 
identified: 

• the elastic parameters: the Young modulus E and the Poisson ratio v; 

• the continuum damage parameters : the limit stress dy and the hardening parameter 

• the discrete damage parameters : the limit normal traction dy, and the softening 
parameter j8 . 

The limits of realistic values for each parameter are shown in the Table [3l Note that 
in our identification methodology we do not suppose to have an expert capable of giving 
the initial estimate of material parameters values, as in e.g. ifTOl I23ll20ll . Therefore, the 
bounds on model parameters were kept rather wide. 

3 Tensile test 

The simplest possibility to identify material parameters for a particular concrete is to per- 
form a uniaxial tensile test. In this case, the stress and strain fields within the speci- 
men would remain homogeneous until the final localized failure phase, and the behav- 
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ior on the structural level is very close to the response of a material point. The load- 
displacement diagram consist of three easily recognizable parts, as shown in Figure [1^: 
the first one corresponding to the elastic response of the material, the second one describ- 
ing the hardening and the third part the softening regime. The calibration of model pa- 
rameters can be carried out to follow the same pattern: first. Young's modulus E and Pois- 
son's ratio v are determined from the elastic part of the load-displacement diagram, fol- 
lowed by the limit stress Of and the hardening parameter K identification from the part 
of the diagram corresponding to the hardening regime and, finally, the limit normal trac- 
tion df and the softening parameter j8 are estimated from the softening branch. Note 
that for Poisson's ratio identification, one additional local measurement is needed to 
complement the structural load-displacement curve, namely the measurement of lateral 
contraction of the specimen, see Figure [1}5. 




1 ' ^ ' ^ ' ^ ' 1 -0 02 1 • ^ • ^ • ^ • I 

0.1 0.2 0.3 0.4 0.1 0.2 0.3 0.4 

Prescribed displacement (mm) Prescribed displacement (mm) 

(a) (b) 

Figure 1 : Tensile loading test: (a) Load-deflection diagram (b) Evolution of lateral con- 
traction. 

Although this kind of calibration procedure is robust and accurate [fT5l|, the experi- 
ment dealing with a simple tension test is rather difficult, if not impossible, to perform 
in a well-reproducible way. For that reason, we turn to study the possibility of parameter 
estimates by using three-point bending tests, which is much simpler to practically perform 
in a laboratory. 

4 Three-point bending test 

In the case of a three-point bending test the global response of a specimen represented 
by the load-deflection {L-u) diagram for the structure cannot be simply related to three- 
stage material response with elastic, hardening and softening part (see Figure [2k)- Nev- 
ertheless, we assume that it will be still possible to employ the three-stage approach de- 
veloped for the uniaxial tensile experiment. Similarly to the previous case, the solution 
process will be divided into the optimization of elastic, hardening and softening parame- 
ters in the sequential way. Each step is described in detail in the rest of this section. 

Due to lack of experimental data, a reference simulation with parameters shown in 
Table |4] will be used to provide the target data. These are the same values as considered 
for simulation presented in Part I of this paper yj. 
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Figure 2: Three-point bending test: (a) Load-deflection diagram (b) Evolution of expan- 
sion of specimen. 



E 

V 

K 



38.0 GPa 
0.1 

2.2 MPa 
1000.0 MPa 
2.35 MPa 
23.5 MPa/mm 



Table 4: Parameter's values for reference simulation. 

4.1 Identification of elastic parameters 

In the elastic range, Young's modulus and Poisson's ratio are determined using a short 
simulations describing only the elastic response of a specimen. Similarly to the uniaxial 
tensile test, the elastic behavior is represented by the linear part of load-deflection dia- 
gram. To identify both elastic parameters this information needs to be supplemented with 
an additional measurement. In particular, we propose to include the specimen expansion 
Al defined as the relative horizontal displacements between the left and the right edge 
of the specimen (as indicated by arrows in the Figure [3]), or in other words AZ = V2 — vi . 
The reference expansion-deflection curve is shown in Figure [2l). 




890 



20 



890 



"2 



Figure 3: Displacements measured to evaluate the expansion A/ = V2 — vi of the specimen. 
The objective function Fi applicable for the determination of elastic parameters can 
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be defined as follows: 

Fx = {Lref{u) -L{u)f-wx + {Mref{u) - M{u)fw2 \ u = 0.01mm (2) 

The load L and the expansion A/ are quantities depending not only on displacement u, 
but also on the values of material parameters. In particular, at the beginning of the load- 
ing regime (where u = 0.01mm), the important parameters are only Young's modulus 
and Poisson's ration, because the other parameters are not yet activated. The quanti- 
ties with index ref correspond to the values taken from the reference diagram. The cor- 
responding value of weights wi and W2 were calculated using 30 random simulations 
to normalize the average value of each of summation terms in Q to one. It is worth 
noting that all the quantities in the objective function are evaluated for the prescribed de- 
flection u = 0.01 mm, which allows the simulation to be stopped after reaching this value. 
Therefore, the first optimization stage is computationally very efficient. 

For the sake of illustration, the shape of objective function Fi is shown in the Fig- 
ures and |4j5. As shown in this figure, the objective function remains rather wiggly 
in the neighborhood of the optimal value. 




Figure 4: Objective function Fi : (a) Whole domain (b) Detail close to optimal value. 



4.2 Identification of hardening parameters 

Once we have successfully determined Young's modulus and Poisson's ratio, we can con- 
tinue towards the estimate of the elastic limit stress dy (representing a threshold of elastic 
behavior) and the hardening parameter K. In the spirit of the uniaxial tensile test the 
limit stress will be related to the limit displacement Uf at the end of the linear part of the 
load-deflection diagram. The hardening parameter K will then govern the slope of the 
diagram in the hardening regime. In our particular case, the slope s is approximated as 
a secant determined from two distinct points, see Figure [51 There are two contradictory 
requirements for that choice: first, the points should not be too close to My to ensure that 
numerical errors due to pseudo-time discretization do not exceed the impact of K param- 
eter; second, they should be close enough to ay^ to ensure that the specimen does not enter 
the softening regime and df will not be reached. If the second requirement is fulfilled the 
corresponding objective function depends only on values of the elastic limit stress dy and 
the hardening parameter K, because Young's modulus and Poisson's ratio are fixed on the 
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Prescribed displacement (mm) 



Figure 5: Measurements for objective function F2 definition. 



optimal values determined during the previous optimization stage. The particular choice 
objectie function adopted in this work is 

s = (L(m/ + 0.01mm) -L(m/ + 0.005mm)) /0.005mm (3) 

leading to the objective function in form 

2 2 
Pl = {Uf,ref- Uf) W3 + {Sref - s) W4. (4) 

To keep this optimization step efficient, the simulations should again be restricted 
to a limited loading range, where the limit displacement can be related to the value 
of Uf,ref from the reference diagram. Note that during optimization process, there is 
no guarantee that df will be exceeded when subjecting the specimen to the limit displace- 
ment. Such a solution is penalized by setting the objective function value to 10 x A^, where 
N — 2 is the problem dimension, see Figure [6^. Moreover, as documented by Figure [6b, 
the objective function is now substantially noisier than for the elastic case and hence more 
difficult to optimize. 

4.3 Identification of softening parameters 

The last stage of identification involves the discrete damage parameters: the limit normal 
traction df and the softening parameter j8 . Variable df represents a limit of hardening 
of material and the appearance of a macroscopic crack. Determination of displacement wy^ 
corresponding to this event, however, is rather delicate. The most straightforward method 
is based on the comparison of the reference curve with a simulation for a very high value 
of Of. The point where these two curves start to deviate then defines the wanted value 
of Uf, see Figure [8^. A more reliable possibility could be based on a local optical mea- 
surement in the vicinity of notch H or acoustic emission techniques [|3]|. In our computa- 
tions we consider local measurements of notch upper comers displacements V3 and V4, see 
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Figure|7l As demonstrated by graph[8b, the 'local' value of My corresponds to the 'global' 
quantity rather well. 




890 



Figure 7: Displacement measured to express crack opening defined as V4 — V3. 

After the specimen enters the softening regime, the shape of both local and global 
curves is fully governed by the softening parameter j8 . Therefore, its value can be fit- 
ted on the basis of the load corresponding to the deflection for which the softening is 
sufficiently active. This leads to the last objective function in the form 

Ft, = {uf^ref — UfYw^-\- {Lref{u)—L{uff'W(, \ u = 0.\5mm. (5) 

Since all the other parameters are already fixed on the values determined during the pre- 
vious optimization steps, this objective function depends again only on two parameters: 
the limit normal traciton df and the softening parameter j8 . 

By analogy to the procedure described in the previous section, the penalization is 
applied to the cases where the specimen does not reach the softening regime until the end 
of simulations. This effect is well-visible from the surface of the objective function shown 
in Figure m 



5 Optimization method 

The algorithm used for the optimization of objective functions Fi, F2 and F3 is based 
on an efficient combination of an artificial neural network, namely the radial basis func- 
tion network (RBFN), and an evolutionary algorithm GRADE extended by the niching 
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(a) (b) 

Figure 8: Comparison of diagrams with and without the spring of crack (a) Load- 
deflection diagram (b) Evolution of difference between chosen local displacements during 
the loading. 




Figure 9: Objective function Fy. (a) Whole domain (b) Detail close to optimal value. 

strategy CERAF. The principle of the algorithm is the replacement of an objective func- 
tion by a neural network approximation and its subsequent optimization by an evolution- 
ary algorithm. The approximation is constructed on a basis of interpolation of several 
points, where the values of objective function were calculated exactly. The approxima- 
tion is adaptively improved by new neurons (points), provided e.g. by optima located 
on the previous approximation. 

The main advantage of this methodology is an inexpensive evaluation of the approxi- 
mation, which is repeatedly used during stochastic optimization process. Computationally 
expensive objective function is evaluated only when new neurons are added to the neural 
network. 

5.1 Evolutionary algorithm GRADE extended by strategy CERAF 

Evolutionary algorithms nowadays belong to the most popular optimization tools. Unlike 
the traditional gradient optimization methods, evolutionary algorithms operate on a pop- 
ulation (a set of possible solutions), applying 'genetic' operators (cross-over, mutation 
and selection). The principles of evolutionary algorithm were first proposed by Hol- 
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land f9l . Ever since, the evolutionary algorithms have reached wide application domain, 
see e.g. books by Goldberg [6J and Michalewicz [21 J for an extensive review. In this 
work, we employ a GRADE algorithm (GRadient-based Atavistic Differential Evolution) 
developed in [ITTIl . The algorithmic scheme of the method is briefly summarized below. 

The first step is to generate a starting generation of chromosomes by choosing ran- 
dom values of all state variables. The size of the population is in our case set to ten 
times the number of optimized variables. After the initiation stage the following steps are 
repeated until a termination condition is reached: 

1 . MUT AT I ON, which creates a new solution x using the operation 



where y is a solution from actual population, z is a randomly created solution and k 
is a real number from given bounds. Ten percents of solutions in a new population 
are created by this operator. 

2. CROSS -OVER, which creates a new solution x according to 



where opt(y, z) denotes a better solution from two randomly chosen vectors y and z, 
c is a real random number with a uniform distribution on the interval [0;2] and s 
changes the direction of the descent (y — z) favor better solution from vectors y 
and z. This operator creates 90% of individuals in a new population. 

3. EVALUATION computes the objective function value for each new solution. 

4. TOURNAMENT SELECTION, where the worst individual from two randomly se- 
lected solutions is deleted. This operator is repeated until the number of solution is 
the same as at the beginning of the cycle. 

The basic version of GRADE is complemented with CERAF strategy flU in order to 
increase the algorithm robustness when dealing with multi-modal problems. The CERAF 
strategy works on the principle of multi-start. If the best value found by a stochastic algo- 
rithm does not change for more than a prescribed value during a certain number of genera- 
tions, CERAF store the optimal point in memory and surrounds it by a 'radioactive" zone 
defined as n-dimensional ellipsoid with each diameter equal to 50% of the size of searched 
domain. If a new solution is located inside the 'radioactive" zone, it is replaced with a 
random one outside the zone. If this solution was created by the CROSS -OVER operator, 
the size of the 'radioactive" zone is decreased by 0.5% of the actual size. During the com- 
putations presented herein, the new local extreme was marked when the best value found 
had not changed for more than the stopping precision (defined as a property of objective 
function) during 100 generations. An interested reader is referred to [Sj for more detailed 
description of CERAF. 



x = y + %-z) , 



(6) 



x = opt(y,z)-Fc5(y-z) , 



(7) 
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5.2 Radial Basis Function Network 



Artificial neural networks (ANNs) were initially developed to simulate the processes 
in a human brain and later generalized for many problems like pattern recognition, dif- 
ferent approximations and predictions, control of systems, etc [i27l l26ll . In this work, 
they will be used 'only" as general approximation tools. The particular implementation is 
based on the idea of radial basis function networks as proposed e.g. in [|22l[T3ll . 

Neural network replace the objective function F{x) by approximation F(x) defined 
as a sum of radial basis functions multiplied by synaptic weights, see Figure [TOl 



c i,b i(x) 




C n,b n(x) 



Figure 10: An approximation using RBFN 



or 

N 

F{x)^F{x) = Y,biix)w,, (8) 

where x is a vector of unknowns, &, (x) is a basis function associated with z-th neuron, w, 
is a weight of the z-th neuron and is the total number of neurons creating the network. 
The basis function bi has the following form 

Z7,(x) = e-ll''-<^'ll'/'-, (9) 

where C/ is a vector of the center coordinates for the z-th basis function and r is a normal- 
izing factor set to 

where dmax is the maximal distance within the domain and D is the number of dimensions. 
Synaptic weights are computed from the condition 

F(c,-)=F(c,-), (11) 

imposing the equality of the approximation and objective function values Yi in all neurons. 
This leads to a minimization problem in the form: 

A' A' 

min£[()^- £^;(c,)w,)2 + A,w2] . (12) 
i=i i=\ 



12 



F 


Stopping precision 


Successful runs 


Maximal number 


Average number 








ot s evaluation 


ot s evaluation 


Fi 


10-3 


100 


32 


16 


Fi 


10-^ 


94 


140 


29 


Fi 


10-3 


80 


140 


47 


F3 


10-^ 


92 


140 


37 


F3 


3.10-3 


76 


143 


47 



Table 5: Summary of reliability study. 



where A,- is a regularization factor set to 10^^. The solution of (fT2l) leads to a system 
of linear equations determining the values of synaptic weights w. 

At this point, the RBFN approximation of the objective function is created and the 
above-mentioned evolutionary algorithm GRADE is used to locate the 'approximate' 
global optima. In the next step, to improve the quality of approximation, three new points 
are added into the neural network: the optimum of previous approximation, a random 
point and another point in the descent direction defined by optima found in two previous 
steps. A more detailed description of the method is available in [14]. 

In the current contest, one additional modification is introduced to deal with the penal- 
ization introduced in Sections 14^21 and |431 If the point to be added to the neural network is 
penalized, it is accepted only in the case, when it was the optimum found by evolutionary 
algorithm. The rationale behind this modification is to avoid the penalty-imposed discon- 
tinuities in the objective function, which are impossible to capture by the smooth RBFN 
approximation. 

6 Identification procedure verification 

When assessing the reliability of the proposed identification procedure a special care must 
be given to stochastic nature of the optimization algorithm. The analysis presented herein 
is based on the statistics of 100 independent optimization processes, executed for each 
objective function. The termination criteria were set to: 

• the number of objective function evaluations exceeded 155; 

• the value of objective function smaller than a stopping precision was found. 

A particular optimization process was marked as 'successful', when the latter termination 
condition was met. Note that since the reference simulation instead of experimental data 
are used, the optimum for every objective function is equal to zero. Results of the per- 
formance study are summarized in Table [5] showing the success rate related to a stopping 
precision together with the maximum and average number of function evaluations. 

Moreover, the different values of stopping precision allow us to investigate the re- 
lation between the accuracy of identified parameters and the tolerance of the objective 
function value. Table[6] shows a concrete outcome of such an analysis, where the maximal 
and average errors are calculated relatively to the size of the interval given by limit values 
for each parameter. 
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Parameter 


Stopping precision on F 


Average error [%] 


Maximal error [%] 


E 


10 ^ 


0.41 


1.23 


V 


1 rv— S 
10 ^ 


0.16 


2.20 


^/ 


10 


0.87 


^ CO 

2.58 


A 




U. /o 


z.4y 




10"^ 


0.30 


0.59 


K 


10-3 


0.49 


1.54 




10-^ 


0.47 


1.32 




10-2 


2.34 


12.21 




3 X 10-3 


0.33 


0.67 




3 X 10-3 


0.26 


2.68 



Table 6: Influence of stopping precision on accuracy of identified parameters. 



The results show that the maximal error for the elastic parameters E and v is less 
than 3%, which is sufficient from the practical point of view. This is also documented by 
Figure [m where no deviation of the reference and resulting curves is visible in the elastic 
range. For the hardening stage (parameters dy and K) a similar precision is unfortunately 
not sufficient as shown by Figure [TTb. Increasing the stopping precision to 10^^ reduces 
the error on parameters roughly by 50%, which is sufficient to achieve almost perfect fit 
of the reference curve. 

Finally, a similar conclusion holds for the parameters 6f and j8 describing the soften- 
ing part of the experiment. The stopping precision equal to 10 ^ is too coarse to achieve 
sufficient precision on parameters and needs to be reduced to 3 x 10 3. The effect of in- 
creased accuracy is then well visible in Figure [TTb . This step completes the verification 
of the algorithm. 
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Figure 11: Comparison of load-deflection diagrams: (a) Hardening parameters (b) Soft- 
ening parameters. 



One more analysis has been performed to estimate the importance of the mesh refine- 
ment for the material parameter identification. Namely, we have repeated the analysis 
on the refined mesh, by using both the exact parameter values and the estimates provided 
by the procedure proposed herein. The results are shown in Figure [l2l where we plot 
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the force-displacement diagram for coarse mesh and fine mesh based computations. We 
can see that the mesh refinement reduces the errors produced with the parameter estima- 
tion. The latter can be combined with the adaptive mesh procedure, which is currently 
being explored [|2]|. 
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Figure 12: Comparison of load-deflection diagrams: computations with both coarse and 
fine mesh, for exact and predicted parameter values. 



7 Conclusion 

We have proposed a sound identification procedure for material parameters of the constitu- 
tive model for representing the localized failure of massive structures. The most pertinent 
conclusions can be stated as follows: 

i) The sequential identification approach employed for the uniaxial tensile test can be 
extended to the three-point bending test. The resulting algorithm is very straightforward 
and has a clear link with the structure of the constitutive model. Moreover, each of three 
stages uses only a part of the test simulation, which leads to substantial computational 
time savings. 

ii) Due to the physical insight into the model, it was possible to construct simple 
objective functions Fi, F2 and F3 with a high sensitivity to the relevant parameters. This 
led to non-smooth and non-convex objective functions, which were optimized by robust 
soft-computing methods. 

iii) The proposed identification procedure was verified on 100 independent optimiza- 
tion processes executed for each objective function. In the worst case, the reliability 
of the algorithm is 76% due to very small number of objective functions calls set in the 
termination condition. From our experience with evolutionary algorithms Q, such a re- 
sult is rather satisfactory. 

iv) As the result of a sequential character of the identification procedure, the errors 
in identified parameters accumulate. Therefore, the values need to be determined with 
higher accuracy then usually required in applications (i.e. 5%) and achievable by neural 
network-based inverse analysis ||T6l . 
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v) The major difficulty of the proposed methods is to properly identify the three stages 
of structural behavior. From the point of view of method verification, where the reference 
load-deflection diagram is not noisy, the problem was successfully resolved. To fully 
accept the procedure, however, the experimental validation of the method appears to be 
necessary. 
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